%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Solveur primal assemblé directe de poutre 1D
%	solveur directe avec une méthode primale du probleme de treaction 
%	d'une poutre
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function [sol,condiSp] = Solveur_primal_assemb (nombreElements,nombreDomaines, F, donnees)


keq = donnees.E*donnees.S*nombreElements*nombreDomaines/donnees.L;

% creation des matrices blocs
Kii = keq*2*eye(nombreElements-1);
Kiidess = -keq*eye(nombreElements-2);
Kii(2:end,1:end-1) = Kii(2:end,1:end-1) + Kiidess;
Kii(1:end-1,2:end) = Kii(1:end-1,2:end) + Kiidess;

Kbb = keq*eye(2);

Kbi = zeros(2,nombreElements-1);
Kbi(end,end)= -keq;
Kbi(1,1)= -keq;

Vtmp = Kii\Kbi'; % resolution d'un systeme pour eviter d'inverse
Ktmp = Kbi*Vtmp;

% creation du complement de Schur primal elementaire
Spelem = Kbb-Ktmp; 

bp=[zeros(nombreDomaines,1);F];

% assemblage du complement de Schur primal
Sp = Assemblage(Spelem, nombreDomaines);

% calcul de son conditionnement
condiSp = cond(Sp);

% pris en compte des C.L. par multiplicateurs de Lagrange
C = [1,zeros(1,nombreDomaines)];
SpLag=[Sp,C';C,0];
bpLag=[bp;0];

% resoltution du systeme 
Uglob = SpLag\bpLag;


% reconstruction de la solution 
sol= zeros(nombreElements*nombreDomaines+1,1);

for i = 1:nombreDomaines
    sol((i-1)*nombreElements+1)=Uglob(i);
    Uinter = (-Kii\Kbi')*Uglob(i:i+1);
    
    sol((i-1)*nombreElements+2:i*nombreElements)=Uinter;
    
end
sol(end) = Uglob(end-1);

sol;
plot(sol);

end
